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Abstract 


Previous simulations studying winds only focus on the line force due to photons from central active galactic nuclei. 
What properties of the winds will be when including the re-radiation force due to the scattered and reprocessed 
photons (i.e., the re-radiation effect)? We perform simulations to study the large-scale dynamics of accretion disk 
winds driven by radiation line force and re-radiation force. For the fiducial run, we find that the re-radiation force 
drives stronger outflows during the early stages. When the flows get into the steadiness, the UV radiation due to 
spectral lines dominates total radiation and the re-radiation effect could be negligible. The opening angle of winds 
narrows as the initial gas density increases. The larger the gas density is, the stronger the re-radiation effect will be. 
For Mgy = 10°M 5, € = 0.3, the outflows do become much stronger with the re-radiation effect and the winds still 
cannot escape from gravitational potential. We find that the detection probability of ultra-fast outflows and the 
properties of the winds are both consistent with the observations. 
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1. Introduction 


Observations evidences of low luminosity active galactic 
nuclei (AGNs) (Crenshaw et al. 2003; Tombesi et al. 
2010, 2014) and black hole X-ray binaries (Neilsen & 
Homan 2012; Homan et al. 2016) show that winds are 
common phenomenon. It is believed that winds play a 
momentous role in AGNs feedback. Both winds and AGNs 
radiation can interact with their host galaxies (Zhang et al. 
2021; Lin et al. 2022). The radiation from AGNs can heat the 
interstellar medium (ISM) and affect the properties of ISM 
(Ciotti et al. 2009; Yuan & Narayan 2014). The winds can 
interact with ISM effectively and result in their co-evolution 
(Ostriker et al. 2010; Weinberger et al. 2017; Yuan et al. 2018; 
Lin et al. 2022). 

Accretion flows of a black hole can be divided into hot 
accretion flows and cold accretion flows. Hot accretion flows 
are highly hot and optically thin, and found at lower mass 
accretion rates systems. Winds from hot accretion flows have 
been directly detected in the low-luminosity AGNs (Almeida 
et al. 2018; Ma et al. 2019; Shi et al. 2021) and black hole 
X-ray binaries (Homan et al. 2016; Munoz-Darias & Jimenez- 
Ibarra et al. 2019). Research has found that the winds from hot 
accretion flows are emitted and accelerated under the combined 
action of centrifugal force, gas pressure gradient and magnetic 
pressure gradient (Yuan et al. 2012; Yuan et al. 2015; Narayan 
et al. 2012; Bu & Gan 2018). 

For the case of the cold accretion mode, it is associated with 
a relatively high mass accretion rate and comprised of cool 
optically thick gas. Observational evidences of the winds from 
cold thin disks have been reported for luminous AGNs 


(Crenshaw et al. 2003; Tombesi et al. 2010; Tombesi et al. 
2014; Liu et al. 2013; Gofford et al. 2015; King & 
Pounds 2015; He et al. 2019) and the soft state of black hole 
X-ray binaries (Neilsen & Homan 2012; Diaz Trigo & 
Boirin 2016; Homan et al. 2016; You et al. 2016). 

Theoretically, there are three important mechanisms, namely 
the thermal driving, magnetic driving, and radiation line force 
driving, which are believed to accelerate the disk winds for a 
cold standard thin disk. The thermal driving mechanism can 
commendably understand the winds observed in some low- 
mass X-ray binaries (Begelman et al. 1983; Woods et al. 1996; 
Higginbottom et al. 2017). The magnetic driving mechanism 
can effectively explain the winds in many astrophysical 
environments (Blandford & Payne 1982; Cao 2014). The 
radiation line force driving is a compelling mechanism to 
interpret the ionization state of the winds and understand fast 
outflows and ultra-fast outflows (UFOs) (Proga et al. 2000; 
Nomura et al. 2016; Nomura & Ohsuga 2017; Yang et al. 
2021). 

For the radiation line force driving mechanism, if the gas is 
moderately ionized and can interact with the UV continuum 
through many UV line transitions, the radiation pressure force 
due to spectral lines (i.e., line force) would accelerate the matter 
and launch the cold winds. The line force can be 1000 times 
larger than the radiation pressure force due to Thomson 
scattering (Stevens & Kallman 1990; Liu et al. 2013), and leads 
to high-velocity disk winds. Hydrodynamic simulations of winds 
driven by the line force have been implemented (Proga et al. 
2000; Proga & Kallman 2004; Nomura et al. 2016; Nomura & 
Ohsuga 2017; Nomura et al. 2020; Yang et al. 2021). It is found 
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that the line force can effectively drive high speed winds from 
the inner region of a cold thin accretion disk and the UFOs 
(Tombesi et al. 2010) are commendably explained by line-driven 
disk wind model. Yang (2021b) perform two-dimensional 
magnetohydrodynamic (MHD) simulations to understand the 
formation and acceleration of UFOs in radio-loud AGNSs. 
Higher-velocity winds with a broader opening angle are 
produced in their MHD models and UFOs in radio-loud AGNs 
may be driven by the combination of line force and magnetic 
field. Ajello et al. (2021) study gamma-rays from fast black hole 
winds. They find that these outflows transfer 0.04 % of their 
mechanical power to 4-rays and the y-ray emission is from 
galaxies with UFOs. 

In order to know the propagation and the properties of the 
winds at large scale (around parsec scale), Cui et al. (2020); Cui 
& Yuan (2020b) study the large-scale dynamics of the winds 
through simulation and analytical methods. They find that 
when the winds from a cold thin disk have a negative Bernoulli 
parameter, the winds will only propagate a short distance from 
their simulated injection inner boundary. However, they 
neglect the radiation force in their simulations. Zhu et al. 
(2022) include both radiation force and radiative cooling/ 
heating and perform numerical simulations to study the large 
scale dynamics of winds. They indicate that the properties of 
the winds depend on the mass of black hole and the luminosity 
of accretion disk. 

However, all the above mentioned radiation driven winds 
simulations only consider the radiation force due to electron 
scattering and UV line processes due to the first order radiation 
coming from the central AGNs. They neglect the scattered and 
reprocessed photons. In reality, the locally generated photons 
through scattering, bremsstrahlung and line radiation in the gas 
around parsec scale can also play an important role. The 
radiation force due to the locally generated photons is called the 
re-radiation force and the resulting impact corresponds to the 
re-radiation effect. 

Sim et al. (2010) study multi-dimensional radiative transfer 
simulation to calculate the spectrum of the hydrodynamics 
simulation of the line-driven accretion disk winds. They find 
that the secondary photons ionize the winds materials. 
Higginbottom et al. (2014) perform Monte Carlo simulations 
of the radiative transfer. It is found that scattered photons can 
also ionize gases, resulting in excessively high levels of 
ionization. However, the postprocess approach, which is 
decoupling radiative transfer and HD/MHD calculations, is 
implemented in their works. Liu et al. (2013) investigate the 
effect of re-radiation force on the dynamics of the flows at 
parsec scale under the irradiation from the central AGNs. They 
find that the accretion flows properties change significantly due 
to the re-radiation effect and outflows become stronger 
consequently. Mosallanezhad et al. (2019) improve the 
calculation of radiation force due to local process and 
make comparison between total radial directions forces and 
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re-radiation forces. They find that the mass flux and velocity of 
the outflows are significantly increased and the outflows are 
mainly driven by the re-radiation force. Obviously, the re- 
radiation effect has an important impact on the accretion flows 
at large scale. However, those simulations concentrate only on 
the large scale of the accretion flows, and they do not consider 
the accelerating of the matter in the inner scale due to the line 
force. Substantial radiative accelerations occur within small 
scale and even lead to high-velocity (0.1-0.2c) disk winds 
(Proga et al. 2000; Proga & Kallman 2004; Nomura et al. 2016; 
Nomura & Ohsuga 2017). They will exchange their momentum 
and energy with the gas. When the winds arrive at parsec scale, 
thus, it is interesting to study the dynamics of the winds driven 
by the combination of radiation line force and re-radiation 
force. 

In the paper, we conduct numerical simulations to 
investigate the large scale dynamics of winds driven by the 
line force and the re-radiation force. In our simulations, 
radiation line force, re-radiation force, and radiative cooling/ 
heating are all included. We first simulate small scale winds 
driven the radiation line force in a dynamical range of 
30-1500r, and the time-averaged properties of winds. Then, 
we perform our large scale simulations in a dynamic range of 
1500-10°r,. We concentrate on the large scale dynamic of the 
winds driven by the line force and the re-radiation force. The 
aim of carrying out small scale simulations is to derive internal 
boundary conditions for the large scale simulations. Winds are 
injected into the inner boundary of large scale simulations. The 
properties of the injected winds come from the time-averaged 
value of winds at outer boundary of the small scale simulations. 

The paper is organized as follows: Section 2 describes our 
models and methods; Section 3 introduces our results; 
Section 4 gives discussions, and a summary is presented in 
Section 5. 


2. Model and Method 


We perform axisymmetric two-dimensional radiation hydro- 
dynamic simulations by PLUTO code, which solves systems of 
conservation laws using the finite volume or finite difference 
approach based on Godunov-type schemes (Mignone et al. 
2007; Mignone et al. 2012; Yang et al. 2021; Yang 2021b; Zhu 
et al. 2022). The basic scenario and physical setup of our 
models are mostly the same as one used in Zhu et al. (2022). 


2.1. Basic Equations 


We carry out two-dimensional numerical hydrodynamic 
simulations. The spherical coordinates (r, 0, à) are employed 
and the equations are solved as follows 


OF oF ap 0, (1) 
dt 
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p= -YP ~ pV® + pRa (2) 
d 
zB — -PV -v + p. (3) 
dt\ p 


Here, p, v, P, e, Ç and ® are density, velocity, gas pressure, 
internal energy, the net cooling/heating rate and gravitational 
potential, respectively. F,44 is the radiation force per unit mass 
including Compton scattering, line force, and re-radiation 
force. We apply an equation of state of ideal gas P = (^y — Le 
with ^ — 5/3. In this work, we consider the gravitational 
potential of both the central black hole ®gy= — GMpgy/ 
(r— r,) (Paczyásky & Wiita 1980) and the host galaxy 
DGalaxy = o, Inr + C (Dye et al. 2008; Bu et al. 2016; Cui 
et al. 2020), where Mag, G, Fs, r, C and o, are the mass of 
black hole, the gravitational constant, the Schwarzschild radius, 
the distance from a point to the black hole, the constant, and the 
velocity dispersion, respectively. We adopt the velocity 
dispersion o, — 2 x 10" ems f. 


2.2. Model Setup 


It is assumed that a geometrically thin and optically thick 
disk is around the central black hole for luminous AGNs. The 
thin disk is responsible for the UV emissions and the hot 
corona is responsible for the X-ray emissions. The hot corona 
can irradiate the thin disk, and the radiation intensity from the 
surface of the accretion disk is written as 


3GMguM, [ r 1/2 
meo = Se (S [1 — (2) 


f 21/2 
PES oio: shi (3) | f (4) 
IT rp TD rp 


Here, rp, rą and M, denote the radial position on the disk 
surface measured from the central black hole, the inner edge of 
the disk, and the mass accretion rate of the accretion disk, 
respectively. f; is the ratio of the X-ray luminosity Ly from the 
corona to the disk luminosity Lp (fy = Lx/Lp). We refer to 
Yang et al. (2021) for the detailed calculation of the value of f;. 
The effective temperature of the disk surface is determined by 
Tar = (rIp(rp)/o) ^, with c being the Stefan-Boltzmann 
constant We consider that the radiation from the high- 
temperature (>3000K ) region of the accretion disk contributes 
to the line force. 

The gas above the accretion disk is irradiated by the UV 
photons from the accretion disk surface and the X-ray photons 
from the corona. Our simulations are implemented in a 
spherical coordinate (r, 0, p). We put the surface of the 
accretion disk on the plane of 0 = 7/2. The plane of 0 = 7/2 is 
above the equatorial plane, at the height z, corresponding to the 
disk scale height. For simplicity of calculation, we assume that 
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the disk scale height is a constant of radius (z = 3.ler,, with € 
being the Eddington ratio, Nomura & Ohsuga 2017). 

The line force and the force multiplier M, depend on the 
ionization of the winds and the spectral energy distribution of 
the radiation field. For more details of the calculations, readers 
are referred to Proga et al. (1998) and Proga et al. (2000). 

We assume that initially the gas above the thin disk is in 


hydrostatic equilibrium (p'(r, 0) = p'(r/2)exp(— —— 2) 
and locally isothermal (T(r, 0) = T.g(r sin 0)) in small scale 


simulations. The density of the surface of the accretion disk 
p'(x/2) is set according to the standard thin disks model 
(Shakura & Sunyaev 1973; Nomura & Ohsuga 2017). For the 
large scale simulations, we do not consider an accretion disk at 
the midplane. Our computational domain covers from 1500r, to 
1.5 x 10°r,. We employ non-uniform grids to discretize the 
computational domain. Our resolution is 300 x 160. Initially, 
low density gas po is put in the computational domain. Note, 
various po, € and Mgg are adopted. The details of various 
models are listed in Table 1. We employ axis-of-symmetry 
boundary conditions at 0 = 0° and adopt equator-of-symmetry 
boundary conditions at 0 = 90°. Refer to Zhu et al. (2022) for 
more details on computational domain and boundary 
conditions. 


2.3. Radiative Cooling/heating 


The radiative processes and the interaction between the 
radiation include Compton heating/cooling, X-ray photoioni- 
zation and recombination, bremsstrahlung, and line cooling 
(Blondin 1994; Proga et al. 2000). The net cooling rate ¢ 
denotes as follows 


pC = n? (Gcompton + Gx — Ly), (5) 
where n is the number density, 
n = p/um,, (6) 
Gcompton 18 the rate of Compton heating /cooling, 
Gcompton = 8.9 x 107?6£(Tx — AT), (7) 
Gx is the net rate of X-ray photoionization heating-recombina- 
tion cooling, 
Gx = 1.5 «10747701 = TZ Tg), (8) 
and Lp; is the rate of bremsstrahlung and line cooling, 
Loi = 3.3 X 10-9" + [1.7 x 10-55 
xexp(—1.3 x 105/T) 
xETIT 1/2 + 10-47 6. (9) 


In the above equations, u= 1 denotes the mean molecular 


weight, m, is the proton mass, and Tx= 10°K is the 
characteristic temperature of the X-ray radiation. The parameter 
ó is used to control line cooling. We set 6 = 1, which represents 


optically thin cooling. 
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Table 1 

Summary of Models 
Model Cold Disk Winds Potential Re-Radiation Magn E Po 
numbers (Mo) (g cm ^? ) 
(1) (2) (3) (4) (5) (6) (7) 
1 RHD86 bh d bh no 105 0.6 225 x 10714275? 
2 RHDS86 bh + g d bh--g no 108 0.6 2.25 x 1071472752 
3 RERHDS6 bh + g d bh -- g yes 108 0.6 225 x 10-!47?/r? 
4 RHD86_bh + g_d25 bh--g no 108 0.6 10 ?5 
5 RERHDS6 bh + g d25 bh+g yes 108 0.6 10 25 
6 RHDS6 bh + g d23 bh--g no 108 0.6 107” 
7 RERHDS6 bh + g_d23 bh--g yes 108 0.6 107” 
8 RHD63_bh_d bh no 10° 0.3 225 x 10-!47?/p2 
9 RHD63 bh + g d bh -- g no 10° 0.3 2.25 x 10-!47?/r? 
10 RERHD63_bh + g_d bh+g yes 10° 0.3 2.25 x 107142752 
11 RHD63 bh + g d25 bh 4- g no 10° 0.3 1025 
12 RERHD63 bh + g d25 bh 4- g yes 10° 0.3 10°75 
13 RHD63_bh + g_d23 bh--g no 10° 0.3 107% 
14 RERHD63_bh + g_d23 bh+g yes 10° 0.3 107” 
Model Cold disk winds Mw Py Ew 
numbers (10Leaa/c”) (Leaa/ c) (Leaa) 
(1) (2) (8) (9) (10) 
1 RHD86 bh d 0.11377 0.14548 0.021582 
2 RHD86 bh 4- g d 0.11372 0.14542 0.021574 
3 RERHDS6 bh + g d 0.11372 0.14552 0.021606 
4 RHD86 bh + g d25 0.11805 0.14605 0.021522 
5 RERHDS6 bh + g_d25 0.11889 0.14486 0.021557 
6 RHD86. bh + g d23 0.76886 0.55599 0.061511 
7 RERHDS6 bh + g d23 0.96922 0.78302 0.079 268 
8 RHD63 bh d 1.8611 x 107 1.4586 x 10? 1.4797 x 10? 
9 RHD63 bh + g d 1.8585 x 107 1.3447 x 1078 1.2835 x 10? 
10 RERHD63 bh + g d 2.6390 x 107 3.8471 x 10-8 7.7709 x 1079 
11 RHD63_bh + g_d25 3.2168 x 10? 9.1662 x 10? 2.9331 x 10 !6 
12 RERHD63 bh + g d25 8.3898 x 10-7? 4.0696 x 10? 2.3195 x 107? 
13 RHD63 bh + g d23 3.2534 x 10? 1.5052 x 10? 6.9671 x 10^ !6 
14 RERHD63 bh + g 423 4.6444 x 10-8 7.9759 x 10? 4.5078 x 107" 


Note. Column (1): model numbers; column (2): model names; column (3): the gravitational potential; column (4): the radiation force (whether the re-radiation force is 
included); column (5): the mass of the black hole; column (6); the Eddington ratio (€ = Lp/Leaa); column (7): the initial density of computational domain. 
Column (1): model numbers; column (2): model names; column (8)-(10) are the time-averaged mass outflow fluxes, momentum fluxes and kinetic energy fluxes of the 


winds measured at the outer radial boundary (r = 1.5 x 10°r,), respectively. 


The net cooling rate ¢ depends on the density p, the 
temperature 7, the characteristic temperature of the X-ray 
radiation 7x, and the ionization parameter €. Here, the 
ionization parameter £ is expressed as 


f. Lo 
E 2 
nr 


exp(— Tx), (10) 
where 7x is the optical depth for the X-ray from the corona. For 
the large scale simulations, we calculate the optical depth for 
the X-ray 7x from the corona and the optical depth for the UV 
radiation Tyy from the accretion disk as follows 


m6. 0) — nx-- f oxp, Or (0 
15007 


Suv Sanat f oep(r', 0)dr'. (12) 
15007; 


where To,x, To,uy are the optical depths for the X-ray and the 
UV radiation respectively across the entire domain of the small 
scale simulations, c, is the mass-scattering coefficient for free 
electrons and set to be 0.4 g cm? We set ox—o0, for 
£210'ergems ! and ox=1000, for £< IO'ergems ! 
(Proga et al. 2000; Yang et al. 2021). 


2.4. The Radiation Force 


Liu et al. (2013) and Mosallanezhad et al. (2019) indicate 
that the re-radiation effect can significantly affect the properties 
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Figure 1. Time-averaged mass fluxes (left panel), momentum fluxes (middle panel) and kinetic power (right panel) of winds for models 1 and model 2. The curves of 


model 1 and 2 are overlapped. 


of the accretion flow around parsec scale. Therefore, the re- 
radiation force can play an important role in driving winds 
process when the winds spread to large scale, which launch 
from inside hundreds of Schwarzschild radius. 

The total radiation force includes the radiation force Ff from 
the center AGNs and the re-radiation force F7 due to the 


locally produced photons, 
Fa = Ff + FS. (13) 


Here, the re-radiation force F;° is in the z-direction, and is 
described by 


Fe = f “(Se + rl dz 2, (14) 


where 
(15) 
is the source term and it denotes the first-order scattered 


photons of the radiation from the central AGNs (Liu et al. 
2013). Fp and F are the accretion disk flux and the X-ray flux, 


Se = poe(Fp + Fx), 


respectively. 
The radiation force F from the center AGNSs is expressed as 
2 
Fe = "OX, Sb M) + P pg, (16) 
pc c c 


In the above equation, the first term on the right-hand side is the 
force corresponding to photoionization heating-recombination 
cooling F cun (Liu et al. 2013), the second term denotes the force 
for the UV radiation F", which corresponds to Thomson 
scattering and line force, and the third term denotes the force for 
the X-ray radiation F®*. We note that the simplification of 
radiation force from X-ray photoionization heating and recom- 
bination cooling will underestimate the force and this approx- 
imation is feasible because of the radiation force in the radial 
direction being in general dominated by the line force. 

There is a difference here in the calculation of the accretion 
disk flux Fp: we calculate it by the specific vector-valued 
integral (Proga et al. 1998, 2000), while Mosallanezhad et al. 
(2019) adopt point-source approximation. Although the point- 
source approximation approach is feasible around parsec scale, 


the vector-valued integral is more accurate and reliable in our 
simulations. In the paper, Fp and Fx are obtained as follows 


Fp = f to(rv)naa, (17) 
f: Lp 

Fx = =, 18 

i 4nr? (18) 


where n is the unit vector, and Q is the solid angle subtended 
by the disk and central object. 


3. Results 


We perform simulations with different values of initial gas 
density pọ and accretion disk luminosity Lp = eLeag (Leaa being 
Eddington luminosity) to explore the re-radiation effect of these 
parameters. In the paper, we also study the influence of galactic 
potential on the dynamics of the winds. The model parameters 
and the general results are summarized in Table 1. The 
properties of the winds in Table 1 are obtained by time- 
averaging the values at the outer radial boundary. We give the 
time-averaged values on time interval of 0.57,5,—1.0T,,, where 
To» denotes the orbital time at the radial outer boundary. 


3.1. The Effect of the Host Galaxy Potential 


Model 1 and model 2 show the winds without and with 
galaxy potential for Mgy = 105 M... € — 0.6, respectively. In 
Figure 1, we plot the time-averaged mass fluxes, momentum 
fluxes and kinetic power of the winds for model 1 and model 2. 
The mass outflow flux, momentum flux and kinetic energy flux 
are calculated as follows 


r 890 
My(r) = An? f pmax(v,, 0)sin 0d0, (19) 
0o 
. 890 
By(r) = Agr? f p max(v?, 0)sin 640, (20) 
0o 
i 890 |] 3 . 
Ey(r) = Axr? f go max (oj, 0)sin 8d. Q1) 
0o 


In the two models, the mass fluxes keep almost unchanged 
with radius and the momentum fluxes and kinetic power 
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Figure 2. Left panel: distribution of specific gravitational force with radius (black line, green line and red line correspond to the gravitational force of galaxy, black 
hole and their sum, respectively). Right panel: the r dependence of time-averaged kinetic energy along 6 = 71° (1/2(v/ c}, black line), galaxy gravitational energy 


(cà In r/c?, green line) and black hole gravitational energy (GMau/(r — rc, red line) for model 2. 


increase due to the acceleration by line force. The winds do not 
stop at the outer radial boundary of our simulations and can 
move out of the AGNs to interact with the interstellar medium 
of their host galaxies. The kinetic power of the winds is 
significantly higher than 0.5% L,qq, so the winds can provide 
sufficient feedback to their host galaxies (Di Matteo et al. 2005; 
Hopkins & Elvis 2010; Bu & Yang 2021). It is clear that the 
curves of the mass flux, momentum flux and kinetic power are 
almost same for model 1 and model 2. This indicates that 
including the gravitational potential of host galaxy has little 
impact on the large-scale dynamics of the winds. 

The left panel in Figure 2 shows the r dependence of the 
specific gravitational force. The black hole gravitational force 
dominates the gravitational force for r « 1.1 x 10°r,, while the 
host galaxy gravitational force dominates for r 1.1 x 106r, 
The total gravitational force is same with the gravitational force 
due to the black hole when r< 10?r, The right panel in 
Figure 2 shows time-averaged kinetic energy of winds along 
0 = 7Y? (at this angle, the winds move almost radially), galaxy 
gravitational energy, and black hole gravitational energy. It is 
clear that the kinetic energy of the winds is significantly larger 
than their galaxy gravitational energy and black hole gravita- 
tional energy, especially when r > 10°r,. Therefore, the galaxy 
gravitational energy and the black hole gravitational energy can 
be negligible compared to the kinetic energy. Thus, there is 
little effect on the large-scale dynamics of the winds for 
including the host galaxy potential. 

We also study the effect of the host galaxy potential for 
various accretion disk luminosity. From Table 1, as we can see 
that compared to model 8, in model 9, the mass flux at the outer 
boundary slightly decrease and momentum flux and kinetic 
energy flux reduce to about 90%. This indicates that the host 
galaxy potential works and restrains the outflows for 
Mpu = 10°M., € = 0.3. The reasons are as follows. The kinetic 


energy of the winds is rather small and the winds cannot reach 
the outer boundary for the case of Mgy = 10°Mo, ¢=0.3. 
Moreover, the black hole gravitational energy decreases 
compared to model 1 and model 2. Thus the effect with 
inclusion of the galaxy potential is more obvious. However, 
those fluxes measured at the outer boundary are significantly 
small and the winds cannot escape from gravitational potential. 


3.2. The Fiducial Run 


We choose model 3 as the fiducial model in our paper and 
compare it with model 2. The parameters of the fiducial model 
are Magu—10 Ms, £=0.6, pg(r) = 225 x 10-47? /r?. 
Obviously, the fiducial model involves the host galaxy 
potential and the re-radiation force. The only difference 
between model 3 and model 2 lies whether the re-radiation 
force is included. 

In Figure 3, we plot some snapshots of density (color) and 
velocity (arrows) for model 2 and model 3. The upper panels 
correspond to t=0.17,, and the bottom panels correspond to 
t=0.5Təæ. The right panels are for model 2, while the left 
panels are for model 3. The winds go away in the direction of 
0 — 71?. From t=0.17,, to t= 0.57 p, it can be seen that the 
winds spread to the outer boundary and the radial speed of the 
winds is continuously increased. This means that the accretion 
disk winds driven by the radiation mechanism move outwards 
and even can go to outside 1.5 x 10°r,. Comparing the upper 
right panel of Figure 3 to the upper left panel, we can see that 
the velocity of the winds in model 3 is slightly larger than that 
of the winds in model 2. Hence, the distance of winds 
propagation for model 3 is further. The gas density close to the 
equatorial plane decreases. Comparing the bottom right panel 
to the bottom left panel, we can see that the velocity and 
density of the winds in model 3 are almost the same with the 
velocity and density in model 2. This indicates that the re- 
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Figure 3. Comparison between models 2 and model 3. Top left: snapshot of density (color) and velocity (arrows) at t = 0.175, for model 2. Top right: snapshot of 
density (color) and velocity (arrows) at t = 0.17, for model 3. Bottom left: snapshot of density (color) and velocity (arrows) at t = 0.5T,, for model 2. Bottom right: 
snapshot of density (color) and velocity (arrows) at t = 0.5T,, for model 3. 
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Figure 4. Properties of the outflows from 0 to 0.57,, for model 2 and model 3. Left panel: time dependence of the mass fluxes of outflows measured at 1.5 x 10°r, 
(solid black line for model 2, solid red line for model 3) and the mass flux ratio K (dashed line). Right panel: time dependence of radial velocities (solid black line for 
model 2, solid red line for model 3) and densities (dashed black line for model 2, dashed red line for model 3). The radial velocities and densities are measured at 
r=1.5 x 10°r,, 0 = 71°. The ordinate on the left-hand side refers to the solid line, while the ordinate on the right-hand side refers to the dashed line. 
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radiation force drives stronger outflows at the first and after a 
period of time the re-radiation effect could be almost 
negligible. In order to show this point more clearly, we plot 
the properties of the winds from 0 to 0.5 To» for model 2 and 
model 3 in Figure 4. 

The left panel of Figure 4 shows the mass outflow fluxes of 
winds (solid black line for model 2, solid red line for model 3) 
on the range of 0-0.57,,. It is clear that the mass flux for model 
3 is larger than that for model 2 from 0 to 0.12 Top. After a 
period of time (about 0.1575), the mass flux for model 3 is 
same with that for model 2. The dashed line of the left panel in 
Figure 4 shows the mass flux ratio K (the ratio of the mass flux 
for model 3 to that for model 2). It can be seen that the mass 
flux ratio is larger than 1 from 0 to 0.1675, and almost equals to 
1 from 0.27, to 0.57,,. The right panel of Figure 4 shows time 
dependence of radial velocity (solid black line for model 2, 
solid red line for model 3) and density (dashed black line for 
model 2, dashed red line for model 3) at r= 1.5 x 10°r,, 
0 = 71°. It is seen that the velocity and density of the outflows 
have a sharp increase at 0.137,,. The velocity and density of 
the outflows for model 3 start to keep steady at 0.237, and the 
outflows for model 2 maintain a quasi-steady state at 0.347. 
When the quasi-steady is formed, the properties of the outflows 
are same for model 2 and model 3. As Mosallanezhad et al. 
(2019) found, the outflows rate of the accretion flows at the 
outer boundary is enhanced under the irradiation from the 
central AGNs when we include the re-radiation force. The re- 
radiation force works and plays an important role in the early 
stage. Thus, the mass flux for model 3 is larger than that for 
model 2 from 0 to 0.16 To» and even two orders of magnitude 
larger than the mass flux for model 2 at 0.127,,. However, the 
value of the mass flux is extremely small during this period. 
After the quasi-steady state being formed, the winds are driven 
by the combination of line force and re-radiation force and the 
UV radiation due to spectral lines dominates total radiation in 
the simulations. Therefore, the mass flux, velocity and density 
of model 2 and model 3 are almost the same. To prove this 
point, we plot Figure 5. 

Figure 5 shows the structure of time-averaged radiation 
fluxes in units of coef for model 3. The left columns give the 
time-averaged radiation fluxes (FP™, Feine) in the radial 
direction and the right columns give time-averaged radiation 
fluxes (FF, F°) in the vertical direction. Comparing the upper 
left panel to the middle left panel of this figure, we can see that 
the value of time-averaged radiation flux F^'"* due to the UV 
radiation from the AGNs in the radial direction is much greater 
than the value of radiation flux F^?"* due to photoionization 
heating-recombination cooling and X-ray radiation. To specify 
this statement, we plot the bottom left panel of this figure, 
which shows the distribution of the radiation fluxes with radius. 
It is clear that Fe!" is almost eight orders magnitude larger 
than F°P™ along 0 = 71°. This indicates that radiation flux due 
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to photoionization heating-recombination cooling and X-ray 
radiation can almost be neglected and the UV radiation 
dominates total radiation in the radial direction when the flows 
get into the steadiness. 

Comparing the upper right panel to the middle right panel of 
Figure 5, it is seen that the time-averaged value of radiation 
flux Fire due to the UV radiation from the central AGNs in the 
vertical direction is much greater than the radiation flux F7^ due 
to the local photons from the local radiative processes and 
scattered photons. From the bottom right panel of this figure, 
we can see that F!"° is almost nine orders magnitude larger 
than FX along 6 = 71°. Thus, the radiation flux due to the local 
photons from the local radiative processes and scattered 
photons can nearly be ignored and the re-radiation effect is 
almost negligible for Mgy = 105M.;, ¢ = 0.6 when the quasi- 
steady state has been formed. Therefore, the properties of the 
winds are almost the same. The kinetic power of the winds is 
significantly higher than 0.5% L.qq and the winds can provide 
sufficient feedback to their host galaxies. There are deep red 
areas within 10?r, in the middle panels of Figure 5 and bulges 
within 3 x 107r, (red line) in the bottom panels. The sharp 
increases are due to the line force multiplier M. As Zhu et al. 
(2022) found, the winds are accelerated when they move from 
small to large radii inside 10°r,, and the line force multiplier M 
is negligibly small outside ^«10?r,. We note, compared with the 
work of Mosallanezhad et al. (2019), the re-radiation force 
seems unimportant due to the low initial gas density. 


3.3. Initial Gas Density Dependence 


We study the effect of initial gas density (po) in this section. 
Figure 6 shows the comparison of properties of the winds for 
models 3, 5 and 7. From the top left panel, top right panel to 
bottom left panel, we can see that the densities, in the region 
(0 « 45?), are increased. The densities and the velocities of the 
winds, in the region (0 > 60?), are the same. In addition, the 
velocities of the winds in the region (45? < 0< 60°) are 
decreased. This means that the opening angle of winds is 
getting narrower and narrower. To show this point clearly, we 
plot the bottom right panel of Figure 6. The bottom right panel 
shows angular distribution of the winds densities and velocities 
at outer boundary for models 3, 5 and 7. It can be seen that the 
densities for models 3, 5 and 7 in the region (0 « 45?) nearly 
maintain at the level of about 10 76 g cm ?, 10 7? g cm ? and 
10 ? gem ^, respectively. The densities at outer boundary 
exactly correspond to their initial densities. In the region 
(45? « 0 « 60°), the densities and velocities change with angle 
and finally each tend to be consistent. In the region 
(60° < 0 < 80°), the curves of densities are the same, so are 
the velocities curves. The reasons are as follows. The winds go 
away in the angle of 45? < 0 < 80? for Mgy = 105 M.., € = 0.6. 
So models 3, 5 and 7 in the region (0< 45°) have few 
velocities and the densities are almost equal to their initial 
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Figure 5. Structure of radiation flux in units of coef (coef = 30, GMgu M, / 47 ?cr;) for model 3. Top left: time-averaged radiation flux F&^?"* due to photoionization 
heating-recombination cooling and X-ray radiation in the radial direction. Top right: time-averaged radiation flux F;° due to the local photons from the local radiative 
processes and scattered photons originally from the central AGNs in the vertical direction. Middle left: time-averaged radiation flux F^'""* due to the UV radiation 
from the AGNs in the radial direction. Middle right: time-averaged radiation flux pe due to the UV radiation from the AGNs in the vertical direction. Bottom left: the 
distribution of the radiation fluxes with radius in the radial direction along 0 — 71?. Bottom right: the distribution of the radiation fluxes with radius in the vertical 


direction along 0 = 71°. 


densities. Initially, low density gas for models 3, 5 and 7 is put in 
the computational domain with po(r) = 2.25 x 1072/7? g cm^?, 
po-— 10” g cm ? and po = 107” g em ^, respectively. In the 
region (45? « 0 « 60?), the injected winds have only few power. 


The greater the initial density, the greater the optical depth, and 
the greater the wind resistance. Thus, the opening angle of winds 
narrows with the increase of initial density. In the winds region 
(60? < 0 « 80°), the winds go away at 0 — 71° and most of the 
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Figure 6. Comparison for models 3, 5 and 7. Top left: contours of gas density (color) and velocity (arrows) for model 3. Top right: contours of gas density (color) and 
velocity (arrows) for model 5. Bottom left: contours of gas density (color) and velocity (arrows) for model 7. Bottom right: angular distribution of the winds densities 
and velocities at outer boundary for models 3, 5 and 7 (black lines correspond to model 3, red lines correspond to model 5 and blue lines correspond to model 7), the 
ordinate on the left-hand side refers to the solid line, while the ordinate on the right-hand side refers to the dotted line. All the data in this figure are time-averaged. 


kinetic energy flux is in this region. The injected winds have 
enough kinetic power to overcome the winds resistance. There- 
fore, the densities (velocities) of the winds, in the region 
(60? < 0 < 80°), are the same. 

From Table 1, we can see that the significance of re-radiative 
effect depends on initial gas density. For the case of 
Mgu = 105M., €—0.6, when po— 10 7 gem ?, the mass 
flux at the outer boundary for model 4 is approximately the 
same as model 5. But when po = ig g cm ^, the mass flux at 
the outer boundary increases by a factor of about 1.5 from 
model 6 to model 7. For Mgy= 109M.., € —0.3, when 
po = 10” gem ^, the mass flux at the outer boundary 
increases by a factor of about 2.6 from model 11 to model 
12. But when pp — 10 7? g cm 2, the mass flux at the outer 
boundary increases by a factor of about 14 with the re-radiation 
effect. It is seen that the larger the gas density is, the stronger 
the re-radiation effect will be. Liu et al. (2013) also find that the 
larger density will produce the stronger re-radiation force. The 
reason is that the re-radiation force is proportional to the square 
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of the density while the gravitational force and gradient of gas 
pressure are proportional to the density (refer to Equation (14)). 


3.4. Accretion Disk Luminosity Dependence 


Zhu et al. (2022) find that in the case of Mgy = 105M.., 
€ = 0.3, the winds cannot escape from black hole potential. In 
this paper, we also study large-scale dynamics of the winds 
driven by the line force and the re-radiation force for 
Mpu = 10fMo, € — 0.3. The corresponding models are pre- 
sented in Table 1. The mass flux of the injected winds (namely 
the mass flux measured at the inner boundary) is equal to 
0.028370L.aa/ c^. From Table 1, it can be seen that the mass 
flux at the outer boundary increases by a factor of about 1.5 
from model 9 to model 10, but these mass fluxes measured at 
the outer boundary are much smaller than that measured at the 
inner boundary. This indicates that the outflows do become 
much stronger with the re-radiation effect and the winds still 
cannot escape from gravitational potential. 
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Figure 7. The 0 dependence of time-averaged column density (solid blue line), 
column density of gas with 3 < log£/(erg cm s^!) < 6 (dashed blue line), 
maximum velocity of winds with 3 < log £/(erg cm s^!) < 6 (solid red line) 
at the outer boundary for the fiducial model. The red dotted-dashed line, the 
red dashed line, the blue dotted—dashed line and the blue dashed line show the 


lines of Umax = 0.03c, Umax = 0.3c, Ng = 10? cm? and Ny = 107^ cm ?, 


respectively. The ordinate on the left-hand side shows the velocity. The 
ordinate on the right-hand side shows the column density. 


We can see from Table 1 that the significance of the 
re-radiative effect depends on accretion disk luminosity. For 
example, when Mpy = 10°Mo, € — 0.6, the mass fluxes at the 
outer boundary are almost same for model 2 and model 3. But 
when Mgy = 10°M, € — 0.3, the mass flux increases by a 
factor of about 1.5 from model 9 to model 10 (the initial 
density is same). The winds with the re-radiation force in the 
case of Mgy = 106M.., € —0.3 have a more obvious effect. 
This result is easy to understand. The mass flux in the case of 
Mag = 105M.., € — 0.6 1s relatively large. Because of initially 
low density gas, the increase of the mass flux due to the re- 
radiation effect is negligibly small compared to the mass flux in 
the case of Mgy = 105M.., € — 0.6. Thus, the mass flux for 
model 2 are almost same as that for model 3. For the case of 
Mpu = 10°M 5, € — 0.3, the magnitude of the mass flux is 
comparable to the increase due to the re-radiation effect. Hence, 
it has a more visible result. The kinetic power of the winds is 
much lower than 0.596L.qq, therefore, the winds might not 
affect the evolution of their host galaxies. 


3.5. Comparison to Observations 


Outflows from AGNSs are one of the fundamental mechanisms 
for the interaction between the central supermassive black holes 
and their host galaxies. Observational evidences indicate that 
UFOs absorbers have been found in hard X-ray bands (6-9 keV) 
in ~40 percents of nearby AGNs (Pounds et al. 2003; Tombesi 
et al. 2010, 2011, 2012; Gofford et al. 2013, 2015; Laha et al. 
2021). We define the UFOs with the outflows velocity spanning 
v ~ 0.03-0.3c, the ionization parameter £ ~ 10°-10°erg cms | 
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and the column density 10cm? <Ny <S 10?^cm ?. It is found 
that the UFOs may be driven from accretion disk by line force 
(Nomura et al. 2016; Nomura & Ohsuga 2017; Yang et al. 2021; 
Zhu et al. 2022). In the large-scale simulations, we also find that 
UFOs are present. Figure 7 shows the view angle (0) dependence 
of column density and velocity of the winds. The outflows over 
44°_66° angle range satisfy the properties of UFOs in our 
fiducial model. The detection probability of UFOs is around 
36%, which defines as Qugo/47 (Nomura et al. 2016). 

Gofford et al. (2015) perform a systematic analysis to the 
outflow properties in a sample of 51 Sukaku-observed AGNs and 
consider the properties of the associated winds. They find the 
mass flux of UFOs Muro; cx ED and the kinetic power of UFOs 
Euros x L$’. Nomura & Ohsuga (2017) and Zhu et al. (2022) 
perform two-dimensional simulations to study the dynamics of 
accretion disk winds driven by line force. It is found that the 
power law indexes are almost same as those given by the 
observations. Typically, the UFOs are detected that the mass 
outflows rate is ~1074-10°g s~', and the kinetic luminosity of 
winds is constrained between ~10*°-10“erg s^. In our simula- 
tions, the winds found in this paper are driven by both the line 
force and re-radiation force. The mass outflow rate falls between 
1.58 x 10? and 1.35 x 107g cm !, and the kinetic luminosity 
ranges between 2.72 x 10“ and 9.96 x 10 "erg s !. The proper- 
ties of the outflows are consistent with the observations. 
Therefore, the UFOs are quite probably driven by the combination 
of radiation line force and re-radiation force. 


4. Discussions 


In the paper, we gain time-averaged the simulation results at the 
outer boundary of the small scale simulations and inject the winds 
into the inner boundary of our large scale simulations. By using the 
time-averaged method, the properties of the high ionization winds 
can be erased due to their low density. Therefore, in the inner 
radial boundary of the large scale simulations, we mainly inject the 
low ionization winds. Thus, in the large scale simulations, we 
mostly explore the dynamics of the low ionization winds. It is too 
expensive to perfectly study the large scale dynamics of winds 
because of inner radial boundary located at close to the black hole 
and outer boundary located at beyond parsec scale. 

In our simulations, we set 0x = 1000, for €< 10°erg cms! 
and ox = c, for £ 2 10°ergcms ' (Proga et al. 2000; Nomura 
& Ohsuga 2017; Yang et al. 2021; Zhu et al. 2022), which is 
inherent limitation in the current setup of hydrodynamics for 
line-driven disk winds simulations. The abrupt switch due to 
the step function is arbitrary. For example, the X-rays are 
heavily attenuated below £= lO0?erg cms! from ox — o, to 
ox = 100c,, resulting in very effective shielding. In the case 
UV radiation, cuy — 6c, for the entire range of £. This 
processing method is physically inconsistent with the absorp- 
tion of ultraviolet radiation by wind material through bound 
transitions. Proga & Kallman (2004) also explore the situation 
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of ox = o, for all €. Nomura et al. (2020) compare the results of 
these two methods and find that the mass fluxes of winds differ 
roughly by a factor of 2. However, we note that both of the 
treatments are arbitrary. In future work, we plan to study more 
realistic physical models of X-ray and UV opacities. 

Sim et al. (2010) and Higginbottom et al. (2014) indicate that 
scattered X-ray photons also contribute to gas ionization with the 
postprocessed radiation transfer calculations. However, radiative 
transfer and hydrodynamics are decoupled in those papers. The 
simulation of actual flows requires coupling radiation transfer 
and HD/MHD calculations, which are not realistic due to the 
large computational costs. In this paper, we perform the 
hydrodynamics simulations of line-driven winds with the re- 
radiation effect. We include the radiation force caused by the 
scattered and reprocessed photons in the momentum equation. 
Although the scattered and reprocessed photons are included in 
the paper, our treatments for calculating the winds photoioniza- 
tion structure and the interaction between the radiation are 
oversimplified. The next generation code of HD/MHD 
formalism is urgent on principle, which can handle radiative 
transfer in a more realistic and detailed manner. 

In addition, the magnetic field is always present. The inclusion 
of magnetic field will enhance the outflows and produce winds 
with a broader opening angle (Kuncic & Bicknell 2007; Cui & 
Yuan 2020b; Yang et al. 2021;). The standard magnetorotational 
instability (SMRI) has been directly observed in a Taylor-Couette 
cell with independently rotating and electrically conducting end 
caps and is responsible for causing the angular momentum transfer 
(Balbus & Hawley 1991; Wang et al. 2022a, 2022b). It is 
interesting to study the large-scale dynamics of winds driven by 
combining the radiation and magnetic field in future. 


5. Summary 


We carry out two-dimensional simulations to study the large- 
scale dynamics of accretion disk winds driven by radiation line 
force and re-radiation force. We study the properties of the 
winds with the re-radiation effect around parsec scale. It is 
currently impossible to perform simulations with the inner 
boundary close to the black hole and outer boundary around 
parsec scale. Therefore, we divide the space scale into small 
scale and large scale. For the small scale simulations, the 
computational domain is 30r, < r< 1500r, The winds are 
driven by radiation line force. For our large scale simulations, 
the computational domain covers 1500r,<r<1.5 x 10°r,. 
The locally generated photons through scattering, bremsstrah- 
lung and line radiation in the gas can also play an important 
role. The dynamics of the winds at large scale is driven by the 
combination of radiation line force and re-radiation force. We 
inject winds at the inner radial boundary (1500r,) which are 
time-average obtained by the small scale simulations. 

Including the gravitational potential of host galaxy has little 
impact on the large-scale dynamics of the winds. For 
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Mgu = 108M, € = 0.6, the properties of the winds are the 
same. In the case of Mgy = 10°Mo, € = 0.3, the mass flux at 
the outer boundary slightly decreases and momentum flux and 
kinetic energy flux reduce to about 90% with the galaxy 
potential. However, those fluxes measured at the outer 
boundary are too little and the winds cannot escape from 
gravitational potential. 

Taking the fiducial model as an example, we find that the re- 
radiation force drives stronger outflows at the first period of 
time. As Mosallanezhad et al. (2019) found, the outflow rate of 
the accretion flows at the outer boundary is enhanced. When 
the flows get into the steadiness, however, the radiation flux 
due to photoionization heating-recombination cooling, X-ray 
radiation, and the local photons from the local radiative 
processes and scattered photons can almost be neglected and 
the UV radiation due to spectral lines dominates total radiation. 
Thus, the properties of the winds are almost the same when the 
quasi-steady state has been formed. The kinetic power of the 
winds is significantly higher than 0.596L,qq, so the winds can 
provide sufficient feedback to their host galaxies (Di Matteo 
et al. 2005; Hopkins & Elvis 2010; Bu & Yang 2021). 

We find that initial gas density can affect the properties of the 
winds driven by radiation line force and re-radiation force. In our 
models, low density gas is put in the computational domain with 
po) = 2.25 x 10-7 ^ scm, pop — 10°” gcm ^, and 
po = 10? gem ?. The opening angle of winds narrows with 
the increase of initial density. Moreover, the significance of re- 
radiative effect depends on initial gas density. The larger the gas 
density is, the stronger the re-radiation effect will be. We also 
study the accretion disk luminosity dependence. For the case of 
Mpu —105M., € — 0.3, the magnitude of the mass flux is 
comparable to the increase due to the re-radiation effect. Thus, it 
has a more visible result. The kinetic power of the winds is much 
lower than 0.5%Leqq, therefore, the winds might not affect the 
evolution of their host galaxies. 

We compare our simulations results to observations. The 
observed dependence of kinetic power and mass flux on the 
accretion disk luminosity can be well produced Gofford et al. 
(2015), Nomura & Ohsuga (2017), Zhu et al. (2022). For the 
UFOs found in hard X-ray bands, the detection probability of 
UFOs and the properties of the winds in our fiducial model are 
both consistent with observations. Therefore, the UFOs are 
quite probably driven by the combination of radiation line force 
and re-radiation force. 
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